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ICASE Fluid Mechanics 


Due to increasing research being conducted at ICASE in the field of fluid mechanics, 
future ICASE reports in this area of research will be printed with a green cover. Applied 
and numerical mathematics reports will have the familiar blue cover, while computer science 
reports will have yellow covers. In all other aspects the reports will remain the same; in 
particular, they will continue to be submitted to the appropriate journals or conferences for 
formal publication. 
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ABSTRACT 

An enhanced diffusion-reaction reaction system (DRS) is proposed as a statistical model 
for the evolution of multiple scalars undergoing mixing and reaction in an isotropic turbulence 
field. The DRS model is close enough to the scalar equations in a reacting flow that other 
statistical models of turbulent mixing that decouple the velocity field from scalar mixing and 
reaction (e.g. mapping closure model, assumed-pdf models) cannot distinguish the model 
equations from the original equations. Numerical simulations of DRS are performed for 
three scalars evolving from non-premixed initial conditions. A simple one-step reversible 
reaction is considered. The data from the simulations are used (i) to study the effect of 
chemical conversion on the evolution of scalar statistics, and (ii) to evaluate other models 
(mapping-closure model, assumed multivariate /3-pdf model). 
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tract No. NAS 1-1 9480 while the author was in residence at the Institute for Computer Applications in 
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1 Introduction 


Turbulent combustion is an interplay of several complex physical processes: e.g. turbulent 
transport, scalar mixing, chemical reaction, heat release. Turbulent transport is present 
only in inhomogeneous flows where the gradients of the scalar statistics are non-zero. Scalar 
mixing- which is important in isotropic as well as inhomogeneous flows- comprises two mech- 
anisms: the effect of the velocity field, which is to cascade scalar energy to smaller length 
scales; and molecular diffusion, which typically acts at smaller scales to reduce scalar fluctua- 
tions. Chemical reaction is a local phenomenon which converts scalar species from reactants 
to products. In problems of interest in combustion, chemical conversion is accompanied by 
heat release and, hence, density reduction. If the heat release is large enough, the resulting 
buoyancy forces can modify the velocity field in which the scalars reside. If the heat release 
is small, the modification of the velocity field is negligible. 

In the scalar-statistics equations, the terms representing several of the above processes 
need closure modeling. A turbulent combustion model is really a composite of these closure 
models of the individual processes. Often, when calculations from the composite model are 
in poor agreement with either experimental or direct numerical simulation (DNS) data, it 
is difficult to isolate the underlying cause (the model of a specific physical process) of the 
disagreement. For this reason, it is useful to construct sample test problems where only one 
or few of the physical processes are important (other processes are either negligible or or 
completely absent). Data from the numerical simulations of such a test problem can be used 
to construct and validate the individual models. 

The divide-and-conquer approach has proved especially successful in the understanding 
and modeling of the inert scalar mixing process. Eswaran and Pope (1988), and subsequently 
others, performed DNS of inert scalar fields evolving from non-premixed initial conditions in 
constant-density, isotropic turbulence. In this first step, the effect of velocity field on scalar 
mixing is studied, without the complexities of spatial inhomogeneities or chemical reaction. 
One of the chief findings of DNS studies is that (for the range of Reynolds number inves- 
tigated, R\ 30 - 90) the velocity field only affects the timescale of the scalar pdf evolution. 
In appropriately normalized time, the scalar pdf evolution is independent of the velocity 
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field. This insight gained from the DNS studies spawned reasonably successful modeling 
approaches in which the scalar evolution is considered only in the normalized time, thus 
decoupling the velocity field from the scalar field. The models are: the mapping closure 
model (Kraichnan 1989, Girimaji 1992b), the Johnson-Edgeworth translation (JET) models 
of Madnia et al (1992) and the assumed /5-pdf model (Girimaji 1991a). 

Study of a test case which combines scalar mixing with chemical chemical reaction in 
isotropic turbulence should be a useful second step towards the ultimate objective of un- 
derstanding and modeling turbulent combustion. In this work, we analyze the test-case 
equations and propose a statistical model for the scalar (species and temperature) evolu- 
tion. The model equation set is a diffusion-reaction system (DRS). The model accounts for 
the effect of the velocity field via a time-dependent diffusion enhancement coefficient factor. 
The model equations are close enough to the original turbulent-reaction equations that other 
models which decouple the velocity field and scalar-field mixing cannot distinguish between 
the two sets of equations. Numerical simulation of the model DRS equations is computation- 
ally much less intensive than DNS of the original equations and may^erve nearly as useful a 
purpose in understanding and evaluating turbulent combustion models. Simulations of the 
DRS equations are performed to (i) shed some light on the effect of reaction on scalar field 
evolution, and (ii) evaluate the performance of other models (developed originally for inert 
mixing) against reacting scalar data. 

The organization of the remainder of the paper is as follows. In Section 2, the relevance 
of the (enhanced) DRS to turbulent combustion is investigated in more detail. Numerical 
simulations of DRS equations is performed in in section 3. A brief description of the numerics 
employed in the simulations is also provided. Section 4 contains the results of the simulation 
and a discussion of its implications to turbulent combustion modeling. The paper concludes 
in Section 5 with a summary. 
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2 Enhanced diffusion-reaction system as a model of 
turbulent reactive flows 

Consider the mixing and reaction of n scalars (of mass fractions (f ) a , a = l,-- n) in a 
velocity field u(X, t) evolving under Navier-Stokes. The following assumptions are made: i) 
no body force, ii) no radiative heat transfer, iii) scalar diffusion governed by Fickian law, iv) 
low Mach number (hence, pressure is considered nearly a constant), and v) high Reynolds 
number (hence, the viscous contribution to the energy equation is negligible). Subject to 
the above assumptions, the scalar and temperature ( T ) evolution equations are 


t d<j> a , d<j> a 1 _ d r n d(j> a] t mn 

P [ Qi a \r ] [P^ot ov I a5 


dt 


'dXi 1 dXr a dX, 


and 


,dh T dh T , 

p[— + Ui • 


d r i dT , ^ d r rTj-i __ lo 

dXl^dX^^dX^ a ° a dx) ^ a 


W r 


( 1 ) 


( 2 ) 


dt ' ~'dXi 

In the above equations, D a and w a are respectively the Fickian diffusion coefficient and 
the chemical production rate of species ot. The heat conductivity of the mixture is A. The 
enthaly of formation and thermal enthaly of species a are ft® and ft„. The total thermal 
enthalpy of the mixture, (ft T ), is given by 


h T = £ tah T a . 


(3) 


ot=\ 


The pressure (p), density of the mixture (p) and the mass fractions are related through the 
equation of state which could be the ideal gas law. It is assumed that the temperature and 
the chemical source terms can be determined given the thermal enthalpy and mass fractions. 

Subject to homogeniety of the turbulent velocity field, the joint pdf of scalar composition 
and thermal enthalpy, F(ft r , <£), evolves according to (Pope 1985) 


,,dF(h T J) d 

p(ft ,<p) 


dt 


a = 1 


(4) 


In the above equation the conditional diffusion of temperature (©t) and scalar (0,s«) are 
given by 
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and 


e * - < 6 > 
where, the notation (6|r) denotes the conditional mean of b with respect to r. The pdf source 
term of temperature ( Rj ) and scalar (Rfa) due to reaction are given by 


and 


Rt = - E *>. 

a— l 


( 7 ) 


R<j>a — Wa- 


( 8 ) 


While the source terms due to reaction are closed in terms of temperature and scalar compo- 
sition, the scalar diffusion terms need closure modeling. The modeling of these terms require 
that the effect of the velocity field on the scalar field be known, at least approximately. 

2.1 Lagrangian coordinate analysis 


The effect of the velocity field on the scalar field is examined in Lagrangian frame of reference 
(x, t ) which is cartesian. For the sake of this analysis, the velocity field is assumed to be 
isotropic. The Eulerian coordinate X is now nonstationary, curvilinear and nonorthogonal, 
evolving according to 

^4=U[X(x, (),<]. (9) 

The Eulerian conditional diffusion can now be written in terms of the Lagrangian deriva- 
tives: 


_ /JL\n f^2ii r a\ - / dXa dxb 

* a ^dx} Da dx)\ W, dx t dx a 
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(10) 


Define 

c ^~m- (11) 

The correlation between Ci„ and the scalar derivatives in the Lagrangian coordinates is 
discussed in detail in Girimaji (1992 c) for the case of inert scalar mixing. It is pointed out 
that Cia related to the rate at which the turbulence stretches a material surface attached 
to x. The determinant of C ta , |C|, is likely to be larger than unity if the material surface 
is stretched by turbulence and smaller than unity if the surface is shrunk. On an average 
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turbulence stretches material surfaces leading to \C\ being larger than unity. This results in 
the Eulerian scalar derivatives (following a fluid particle) to be larger than ones in Lagrangian 
coordinates. The length scale of variation of (7, a is likely to be of the order of the characteristic 
length scale of the small scales which cause the strecthing or shrinking. This scale is the 
Kolmogorov length scale, 77 . In the absence of reaction, the length scale of variation of the 
Lagrangian scalar derivatives will depend largely on the initial length scale of the scalar 
field. It was argued in M Girimaji (1992c) that, if the Kolmogorov length scale is of a 
different order of magnitude than the initial length scale of the scalar field then, C{ a and the 
Lagrangian scalar derivatives would be poorly correlated. 

Chemical reaction can cause the Eulerian scalar derivatives to be different from those in 
inert mixing. It is important to evaluate the validity of the poor-correlation simplification 
in the reacting case. In the distributed regime of turbulent combustion (Damkohler number 
smaller than unity), the reaction zone (flame thickness) is larger than the Kolmogorov length 
scale. It is unlikely that chemical reaction can cause the scalar derivatives to be significantly 
different from the inert-mixing case. As a result, as in the case of inert mixing, the Lagrangian 
scalar derivatives would be poorly correlated with Ci a . In the flamelet regime of turbulent 
combustion (Damkohler number greater than unity), the flame thickness is smaller than the 
Kolmogorov length scale. Infact, at the limit of infinite Damkohler number, the reaction zone 
can be infinitesimally small, causing a discontinuity in the scalar fields. There is no chemical 
reaction outside these discontinuities. In the reaction zone, it is not clear if the Lagrangian 
derivative will be poorly correlated with (7, a , whereas, outside that zone (regions of no 
reaction) the poor-correlation simplification is as valid as in the inert-mixing case. Inside 
the reaction zone the chemical source term is dominant and the molecular diffusion term- 
and hence its modeling- is unimportant. With increasing Damkohler number - decreasing 
flame thickness - the volume of the flow field over which the simplification is valid increases 
and the importance of the diffusion term inside the reaction zone diminishes further. In the 
limit of Damkohler number going to infinity, the simplification is valid almost everywhere 
in the flow field except at the inifinitesimally thin flame sheet where the molecular term is 
insignificantly small. As a result, the poor-correlation simplification should yield reasonable 
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results even in the flamelet regime of combustion, especially at large Damkohler number. 
Hence, the expression for the conditional diffusion can be simplified as 

t&J'a $ r p, i\ / ^6 i y ^ f 1\ /in\ 

" ** ~ W Mi lh} pDa ~d^}\ ^WiW^lh} pDa l^\ h (12) 

The quantity is a second-order tensor which, if the turbulent field is isotropic, can 

be written as 


<feit> - < s «w- 

where S(t ) = Hence the model for 0^ a is 


e ta (T,f, = {S(t))^{D^}\h T ,V). 


( 13 ) 


( 14 ) 


Using similar arguments the model for Qj can be derived: 


( 15 ) 


The model. Based on the Lagrangian frame analysis, it is proposed that conditional 
diffusion be modeled by 


and 




e *“ = 


( 16 ) 


( 17 ) 


In the above equation, the superscript m refers to the model variables. 

The effect of the velocity field on the scalar evolution in the Lagrangian reference-frame 
is much smaller than the effect in the Eulerian reference frame. Under the poor- correlation 
simplification, the scalar field evolution in the Lagrangian frame can be considered indepen- 
dent of the velocity field. Hence, the model variables ( h Tm and <f > m ) can be considered to 
evolve according to 

+ ( 18 ) 

and 

dh Tm .... a 8T JL „ 8<h m JL . 

( 19 ) 


dh Tm 8 8T n 8(h m n 

pAt- = + E - E *>.(T)- 


dt 
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Equations (18) and (19) have two important uses. First, they form a statistical model of 
the chemically-reacting flow equations (1) and (2). Second, the model equations are similar 
enough to the original equations that many of the other scalar mixing models - multiscalar 
mapping closure models (Pope 1991, Gao and O’Brien 1991, Girimaji (1993)), multivariate 
/?-pdf model, JET model - cannot distinguish the difference between pdf’s of the two sets of 
instantaneous equations. Hence, for validating many statistical aspects of the other models, 
DNS of equations (18) and (19) can be used rather than the more expensive DNS of turbulent 
reacting flow equations (1) and (2). 

What is (S'(t))? The quantity ( S(t )) is a measure of the mean deformation of the Eule- 
rian coordinate frame caused by the velocity field and can be interpreted as the diffusion 
enhancement factor. It can be expressed in terms of the deformation characteristics of a 
material cube attached to a Lagrangian point x (Girimaji 1992c): 

5(x,<) = 3^WW + ^W + ^IW1. (20) 

where A a (x,t) is the area of surface of the cube which was initially coincident with the a 
axis and V(x,t) is the volume of the cube. The initial values are such that, 

-4i(x,0) = ^ 2 (x, 0) = ^ 3 (x,0) = l/(x,0) = S^XjO) = 1. 

For an incompressible velocity field, the volume F(x, <), which also represents the density 
ratio following a fluid particle, is always unity. The compressibility of the velocity field 
manifests itself on the scalar mixing through V. In an isotropic velocity field, the various 
area magnitudes {A\, A 2 , and A 3 ) are all statistically equivalent. Hence, assuming isotropy 
(and incompressibility) 

{S(t)) = (A 2 {x,t)), (21) 

where A(i) represents the area at any timet of any typical material-surface element (of initial 
unit area) associated with the fluid element. Material element deformation in isotropic tur- 
bulence was studied by Girimaji and Pope (1990) and its implications on (S(t)) is discussed 
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in detail by Girimaji (1992c). It suffices to say here that in isotropic turbulence ( S(t )) grows 
in time nearly exponentially. The logarithmic growth rate is approximately 


p(t) = 


dln(S) _ 08 

dt t„ ' 


( 22 ) 


where r, ( is the Kolmogorov time scale of the turbulence. The diffusion enhancement factor 
( S(t )) reflects the cascading effect of the velocity field on the scalar field. As the scalar 
energy cascades from large scales to small scales, the Eulerian scalar gradient following a 
fluid element increases exponentially. The diffusion enhancement effect of the velocity field 
is more prominent at later times when more of the scalar energy is in the small scales. 

3 Further simplification and simulation 

We want to consider the simplest non-trivial case of multiscalar mixing and combustion 
where reaction affects the evolution of scalar pdf. For this purpose, several previous authors 
(e.g., Pope 1991) have used the following set of equations: 

d 2 <f> a 


d <t>a + Vi ^a_ = ^ 




(23) 


dt ' ~'ox, “ dx,ax , 

The above equation represents constant-density, isothermal chemical reactions between scalars 
of different, but, constant diffusivities. If the velocity field U is isotropic, the enhanced- 
diffusion reaction model corresponding to equation (23) is 

(24) 


■W PtfS . ... 


dt '~ K ~"~ a dx a dx a 
The time-dependent diffusion enhancement factor is calculated from equation (22). 

We consider reversible reactions of the type 

A + B ~ 2P, (25) 

where, scalars A and B are reactants and scalar P is the product. The factor two is included 
in equation (25) to conserve mole fractions. The chemical production rate of the three species 


are 


w A = ~kj<f>A<t>B + h4> c 
W B = -kj<f>A<!>B + h<j>l 

wp = 2kf<f>A<t>B - 2 k b <f> 2 c , 


( 26 ) 
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where, kj and k b are the forward and the backward reaction-rate coefficients. The backward 
coefficient is usually expressed in terms of the forward rate: k b = k e k 5 , where k e is the 
equilibrium constant. 

3.1 Simulations 

The problem of non-premixed, randomly and isotropically distributed scalars undergoing 
enhanced mixing and chemical reaction is considered. Numerical simulations of equation 
(24) are performed in a cubical domain with periodic boundary conditions. A pseudo- 
spectral method is used: the scalar spatial derivatives are sought in spectral space, whereas, 
the reaction-rate term is calculated in physical space. The initial scalar field is specified in 
the manner described in Girimaji (1993). The initial scalar field has a prescribed length 
scale 1 $, and the initial joint pdf, F{<j>A , <^b), is given by 

F{<j>A, <f>B ; 0) = //^(O)^^ — 1)<H</>b) + I^B{0)S(<f>A)H < f ) B — 1) + /*p(0)6(<^4)£(<£b)- (27) 

In the above equation, ha, Pb and I^P represent the means of species A , B and P respectively. 
The mass fraction of P is given by 

4>p — l — <t>A — 4 >b ■ (28) 

Even when the scalar field evolution is governed by the simple enhanced-diffusion reaction 
system (equations 24 and 26), the scalar joint-pdf evolution is affected by several parameters: 
i)r„ (the Reynolds number); ii) 1 4 , (initial length scale of the scalar field); iii) D A , D B 
(Prandtl number and relative diffusion); iv) ^ a , Hb (stoichiometric ratio); v) kf (Damkohler 
number); and, vi) k b (reversibility). In this study, we fix the values of T n (corresponding 
to a Taylor-scale Reynolds number of approximately 60, as reported in Girimaji and Pope 
1990). The simulations are divided into four groups based on the initial joint pdf. Within 
each group several simulations are performed by varying the other parameters. The values 
of the parameters for the various simulations are shown in Table 1. 
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4 Results and implications to turbulent combustion 
modeling 

The results of the simulations are presented in two parts. First, the effect of reaction on 
the evolution of the means, variances, pdf’s and other quantities of interest are studied. 
Secondly, the simulation data are used in its role as the facsimilie of DNS to evaluate two 
models of turbulent combustion: the mapping closure model and the assumed /3-pdf model. 
The effect of chemical reaction on the following statistics of various species are examined: 

1. The mean, variances and correlations. 

2. The evolution of the pdf. 

3. The mean scalar dissipation and chemical production rate of variance. 

4. The evolution of the conditional scalar dissipation. 

Results from the simulations of Group 1 and Group 2 are used in this part of the study. 
Statistics of species A and P only are presented, since <j>B can be completely determined 
knowing 4>a and <^p. In the discussions below, angular brackets imply mean and primes 
denote fluctating part of a random variable. 

Means. The mean evolutions for Groups 1 and 2 are presented in Figures 1 and 2 respec- 
tively. As is to be expected, the means do not change for the inert cases: (<j> a) remains at its 
initial value, and ( <f>p ) is always zero. As for the reacting cases, the mean reactant decreases 
more rapidly for larger kj (reaction rate coefficient) than for smaller kj. In the reversible 
reaction cases (Clb and C2b), the mean values asymptote to equlibrium values given by the 
equation 

PaPb = Vp- (29) 

The above equation in conjunction with equation (28) leads to the following equlibrium 
values: for case Clb, \i A = 0.33, hb = 0.33, [ip = 0.33; and for case C2b, /j a = 0.116, hb — 
0.616, up = 0.267. 


10 



Variances. The evolution of variance in reacting field in given by 


= 2(-4 + 0, (30) 

where, the mean scalar dissipation (c^) and the production/destruction of variance due to 
reaction (e^) are defined as 

d = M'o, <W'« . 

a K dXidXi h 

e; = (w'M. (32) 

The variances of species A for Groups 1 and 2 are shown in Figures 3 and 4. In Group 1, 
reaction appears to cause the variance to decay slowly compared to the inert case. On the 
contrary, in Group 2, the variance decays more rapidly with reaction than without. The effect 
of reaction on Group 2 is larger than that on Group 1. The reasons for these observations 
are explained further below, when the behavior of e A is discussed. Although discernable, 
the magnitude of the modification due to reaction is small enough to be unimportant for 
practical problems. As shown in Figures 6 and 7, e d A ^ e T A for the cases considered resulting 
in only a slight modification of the variance evolution. 

The variance of the product is initially zero. It grows in time initially, due to the chemical 
production term; more rapidly in Group 2 than in Group 1. After attaining a peak value, it 
diminishes due to effect of molecular action. 



Correlation between (j) A and The correlation coefficient between species A and B, 


C AB = 


Wb) 




2 _2 
A a B 


(33) 


is given in Figure 5 for the inert cases C2 and C'4 and their reacting counterparts C2b and 
C4b. At early times, when chemical conversion is small, correlation is not affected much 
by chemical reaction. With time the correlation for the reacting cases deviate from their 
inert counterparts, going to lower magnitudes while still preserving a negative sign. The 
deviation, however, is not too large. 
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Scalar pdf. In order to examine the effect of reaction on the scalar pdf evolution, the pdf 
of scalar A is plotted in Figure 6 for cases Cl and Clb. Figure 7 contains pdf’s of scalar A 
for cases C 2 and C2b. Initially, chemical reaction does not affect the pdf much, for, given the 
segregated initial condition, molecular mixing should occur first before chemical conversion. 
However, at later stages, as chemical conversion becomes prevalent, the pdfs are strongly 
affected. The pdf’s develop positive skewness due to chemical depletion, and in both cases 
tend to ^-functions at zero at long times. 

The pdf of the product, scalar P, is plotted in Figure 8 (Cases Cla and Clb) and Figure 9 
(Cases C2a and C2b). The pdf’s of both groups display bimodal behavior initially. One mode 
of high probability is at zero value of mass fraction, representing parts of the field containing 
unmixed reactants. The second mode of high probability is close to the maximum value of 
mass fraction at that time. With time, the mode at zero disappears, and the other mode 
migrates to higher values retaining its spike-like form, in cases Cla and Clb. The asymptotic 
state of the pdf for Group 1 is a ^-function at unity, representing a complete conversion of 
all reactants to products. In cases C2a and C2b also the mode at zero dissipates with time. 
The pdf migrates to the right and has wider support (larger variance) unlike Cla and Clb. 
The asymptotic form of the pdf in these cases (Group 2) is a ^-function at cf>p = 0.5. 

Chemical production-rate of variance. The chemical production-rate of variance of 
species A is given by 

e A = ( w a4‘a) = ~(4a)(4a4b) ~ (4b)(4a 4a) ~ (4a4a4b)- ( 34 ) 

For Groups 1 and 2, the correlation between A and B is nearly negative unity leading to 

*'a » -*B< ( 35 ) 

which when substituted in equation (34) leads to 

= WaPa) ~ 1(4a) - (4b)](4'a4'a) + (# 4 )- ( 36 ) 

Similar arguments lead to 

eg = ( tv' B <f> l B ) ~ [(4b) ~ (4 a)](4'b4'b) + (4b)- ( 3 ?) 
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In Figures 10 and 11, e A of Groups 1 and 2 are presented. The behavior of e A in the two 
groups can be understood by recalling that for Group 1, (phi a) = ( <t>B ) and for Group 2, 
(4>b) — ( <f>A ) + 0-5, at all times. This yields 

<4 * (€)■ (38) 

for Group 1. As mentioned earlier, the pdf is positively skewed leading to a positive value 
for e A . For Group 2 

e A « -0.5 (phi A ) + (<f>'l). (39) 

For this group, the first term on the right hand side (of equation 39) dominates the skewness 
term leading to a larger negative value of e A . The behavior of the variances of the two groups 
is consistent with the above explanations. 

Figures 12 and 13 contain t T P of groups 1 and 2. Both groups exhibit similar behavior. 
Starting from zero, the value of e P increases to a peak value and then diminishes to near-zero 
values. The peak value is higher for higher kj. The early behavior of t r P of the reversible 
case is identical to its non-reversible counterpart. However, its peak value is not as high, 
and it diminishes much more rapidly to zero. 

Mean scalar dissipation. The mean scalar dissipation of species A of Groups 1 and 2 
are given in Figures 10 and 11. Similarly to its effect on variance, reaction causes the mean 
scalar dissipation to decay more slowly in Group 1 and more rapidly in Group 2. Again, 
the magnitude of modification is not very large. In both groups, the magnitude of the mean 
scalar dissipation is much larger than that of chemical production-rate of scalar variance. 

The mean scalar dissipation of P of Groups 1 and 2 are given in Figures 12 and 13. In 
both the groups, with the formation of products, the mean scalar dissipation increases from 
a zero initial value. After attaining a maximum value, the dissipation decreases. The mean 
scalar dissipation achieves higher values in Group 2 than the corresponding cases in Group 
1. Initially, e P is larger in magnitude than e P , leading to a growth of a P from its initial zero 
value. At latter times, e P is larger leading to a decay in the variance levels. The implication 
is that initially the scalar-P evolution is reaction controlled, after which both reaction and 
diffusion are equally important, and ultimately the evolution is diffusion controlled. 
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Conditional scalar dissipation. The evolution of the conditional scalar dissipation of 
scalar A in cases C2 (inert) and C2b (reacting) is given in Figure 14. Chemical reaction 
has two effects on this quantity. The first is due to the more rapid decline in the maximum 
value of <f>A as a result of reaction. The zero of the conditional dissipation migrates with the 
extremum value (Girimaji 1992b). When mixing is accompanied by reaction, the conditional 
dissipation is non-zero over a smaller range of 4>a > than in the case of inert mixing. This 
effect is negligible at the early stages and more pronounced in the later stages. The second 
difference is the higher value of the conditional dissipation (where it is non-zero) at a given 
value of mass fraction in the reacting case as compared to inert case. Again, this difference 
is large in the latter stages and negligible in the early part. Despite the higher values of 
the conditional scalar dissipation in the reacting case than in the inert case, the mean scalar 
scalar dissipation is indeed lower in the reacting case, as shown in Figure 11. The reason for 
this is the shift in the pdf (Figure 11) due to reaction. In the reacting case, the high values 
of conditional dissipation occur at values of mass fraction of low probability of occurence 
and vice versa. Whereas, in the inert case high probability and high conditional dissipation 
appear to occur at nearly same values of mass fraction. 

The conditional scalar dissipation of the product, scalar P, is given in Figure 15 for cases 
C2a and C2b. The conditional dissipation of the product is very unlike that of the reactant 
for this case. Initially, it is non-zero over a very narrow range of <j>p values and it gradually 
widens with time (or reaction). The peak value normalized by the mean scalar dissipation 
decreases with time, indicating gentler gradients in the product field at later times. The 
behavior of this quantity is further discussed later when the modeling issues are examined. 

So far in this section, the effect of reaction on various quantities of interest was examined. 
It is found that although the mean and pdf of the scalars are strongly affected by reaction, 
other quantities like the variance (and, perhaps, other even moments) are not very different 
from their inert-mixing counterparts. Hence, it would be useful to directly compare and 
evaluate models of inert mixing against the enhanced-diffusion/reaction data. 
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4.1 Models vs. Data 


We attempt to answer three specific questions regarding turbulent combustion modeling 
using other models: 

1. Can the inert mixing multiscalar mapping-closure model (Girimaji 1993) be used for 
reacting flows without modification? 

2. How good is the multivariate assumed /?-pdf model (Girimaji 1991a) for calculating 
reacting flows? 

3. How does the computationally simple assumed-pdf model compare with mapping- 
closure model for multiscalar mixing? 


Multiscalar mapping closure model. The multiscalar mapping closure model for inert 
mixing (Girimaji 1993) is based on the simplification that the conditional scalar diffusion of 
a given scalar is a function of that scalar only: 

e - = K ( wFwF M - (40) 


The conditional scalar diffusion for each scalar is obtained by employing the mapping closure 
procedure for single scalar mixing. Knowing the conditional scalar diffusion of each scalar, 
the joint pdf evolution is solved, and the results are in reasonably good agreement with data. 

Even when mixing is accompanied by reaction, it is suggested in Girimaji (1993) that 
the simplification stated in equation (40) may be valid. In the reacting case, however, the 
mapping closure procedure even for single scalar is not clear. So it will be useful to know, if 
the conditional scalar diffusion implied by the mapping closure procedure for inert mixing is 
adequate for the reacting case also. Conditional scalar diffusion is related to the conditional 
scalar dissipation according to 


G (p a 


1 dx*»W«) 


(41) 


F(<f) a ) d<t>a 

The inert case mapping closure model (with Gaussian reference field) for conditional scalar 
dissipation for initially non-premixed reactants is (Girimaji 1992b) 


X(0.5) 


= exp(-2[erf { 2 <f> - 1 }] ) 


(42) 
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This model is quite good in the early stages of inert mixing and is technically invalid but 
still adequate during the latter stages. The validity of the above model for the reacting case 
is now investigated. 

Shown in Figure 16, are F((f>A)x{4 > A)/x{ < f > A — 0.5) of the data and model for case C2b. 
In the model calculation of the above quantity, pdf F is taken from the data. The model 
agrees with the data very well at early times. At later times, the model is still adequate. 
Recall that even for the case of inert mixing the model is not very good at the final stages. 
Perhaps, an inert-case model which is uniformly valid at all times would be good for reacting 
case at all times too. In any event, for combustion applications, the behavior of the model 
at the early stages, when the unmixedness is high, is more important than at the late stages 
when the scalars are more or less uniformly mixed. Hence, for the reactants a closure model 
for the conditional scalar diffusion obtained by substituting equation (42) in equation (41) 
might be adequate. 

Modeling the conditional diffusion of the product (scalar P) is not as simple. The relative 
importance of chemical conversion and mixing in the scalar pdf evolution (of P) at various 
times can be surmised from Figure 9 where the mean scalar dissipation and the chemical 
production rate of variance are compared. Initially, the evolution of the scalar P is dom- 
inated by chemical conversion. Molecular diffusion does not play a significant role in the 
pdf evolution until much later. Therefore, it is important that the model for the conditional 
diffusion (or dissipation) of the product be accurate at later times; accuracy at early times is 
not as crucial. Comparison of F((f>p)x{<t>p) / xi^P = 0-5) f° r the product is shown Figure 17 
for case C2b. The agreement is very poor in the initial stages. At the later stages, it is much 
better. Whether the model is good enough for engineering applications can be determined 
only from the computation of the pdf evolution using the closure model. Such a computation 
is outside the scope of this work. 

Assumed multivariate /3-pdf model. In the assumed-pdf approach the scalar joint pdf 
is prescribed knowing the first few moments of the scalar field. It is computationally far less 
intensive than the mapping closure model and is ideally suited for engineering computations. 
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The assumed multivariate /5 pdf model for a iV-scalar mixing process is given by (Girimaji 
1991a) 

-*> - - - M- < 43 ) 

- l (Pi) • • • HPn) 

The parameters of the model (/5 1? ■ • • , /5/v) are functions of the mean mass fractions fi a and 
turbulent scalar energy Q: 

i _ C 

(44) 


“ >)• 


In the above equation Q and S' are given by 


AT N 

Q = I> 2 ; s= X/^ 2 - 

ar=l a=l 


(45) 


where cr^ is the variance of mass fraction of scalar a. The above model tested for two-scalar 
inert mixing process with reasonable success, but is yet to be validated for multiscalar mixing 
and reaction. 

In Figure 18, the individual variances and some cross moments of species A and B cal- 
culated using the model are compared against corresponding simulation data for case Clc. 
Similar comparisons are performed for cases C2c, C3b and C4b in Figures 19, 20 and 21, 
respectively. The means and the turbulent scalar energy required as inputs to the model are 
taken from DRS data. The model predicts the variances reasonably well. The magnitude of 
the cross covariance is consistently underpredicted by the model, but the difference is not too 
much. For the two higher order cross moments compared {(^a^'b) anc ^ (^a^b)) ^ ie m °del 
does surprisingly well. For each of the cases considered, the moments behave differently and 
the model is able to capture the behavior well, qualitatively and quantitatively. Overall the 
performance of the model is quite satisfactory, given the simplicity of the model. 


Multivatiate /5-pdf model vs. Mapping closure model. Now the assumed /5-pdf 
model is compared against the more detailed and computationally intensive mapping closure 
model. The two models are equally good for the case of inert non-premixed two scalar 
mixing (Girimaji 1992b). Comparison of the two models for multiscalar mixing has not 
been performed before, and is attempted presently for inert mixing. (As mentioned before, 
computations of multiscalar reacting flows using mapping-closure model are intensive to be 
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attempted here.) For the purpose of computation, simulations C3 and C4 are used. In the 
absence of reaction, Groups 1 and 2 reduce to two-scalar mixing, and hence not useful to 
evaluate multiscalar. In Figure 22, various normalized cross moments (of species A and B) 
obtained from the simulations are plotted as a function of the variance of scalar A for Case 
C3. (The variance decreases monotonically in time and hence can be used in lieu of time 
as the abcissa.) Also shown in the Figure are the cross moments calculated using the two 
models. The moments of the normalized variables are presented: 

q — tz. — it°t 

o a 

For all the moments considered, both the models agree equally well with the data. Figure 
23 presents the same comparison for case C4. In this case, the performance of the mapping 
closure model is clearly superior. These findings are in keeping with the arguments presented 
in Girimaji (1991a) that the assumed /3-pdf model is likely to be more accurate when the 
length-scale and diffusivities of the scalars are similar (case C3) than when they are widely 
disparate (case C4). For case C4, even the mapping closure model is not very accurate. This 
leads to the question, for engineering calculations, is the increase in accuracy achieved using 
the mapping closure model worth the extra effort involved in computing? The answer will 
depend on the application for which the model is being used. In combustion calculations it 
is the unnormalized moments that appear in the scalar moment equations and need closure 
modeling, not the normalized moments. Despite the relatively poor agreement of the nor- 
malized moments in Figure 23, the unnormalized moments calculated using the /?-pdf model 
are quite close to the data as seen in Figure 21 for the corresponding reacting case C4b. 
This appears to suggest that the assumed-pdf approach, although simplistic, may be quite 
adequate for modeling of turbulent combustion in engineering calculations. 

5 Conclusion 

It is shown that in isotropic turbulence the (enhanced) diffusion-reaction system (18 and 
19) can be considered a statistical model of turbulent combustion (equations 1 and 2) with 
the initial field and the diffusion enhancement factor ((5(f))) as inputs. The diffusion- 
enhancement factor can be found if the scalar variance evolution is known (Girimaji 1992c). 
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In the event that scalar variance evolution is not known from any other source, a simple 
model for the diffusion-enhancement factor is provided (equation 22). The enhanced-diffusion 
reaction model is close enough to the turbulent combustion equations that, in normalized 
time, many of the other models (mapping closure models, assumed-pdf models and any 
model that decouples the velocity and scalar fields) cannot distinguish between the present 
model and the original equations. 

Simulations of the diffusion-reaction system is performed in a cubical box (64 3 grid points) 
from non-premixed initial conditions. The reaction considered is of the type 

A + B <-+ P. 

The reaction is constant-density and cold (isothermal). Several simulations of different initial 
conditions, reaction-rate coefficient and degree of reversibility are performed (Table 1). The 
simulations are used to examine the effect of reaction on the evolution of several scalar 
staistics (Figures 1 - 15). The data are also used to evaluate other models; mapping closure 
model and the assumed multivariate /5- pdf model. The observations and inferences from the 
study are the following. 

1. The conditional scalar dissipation implied by the inert mapping closure model is rea- 
sonable even for reacting case but only for the reactants (scalars A and B) as shown in 
(Figure 16). For the product (scalar P), the agreement is quite poor at early times but 
better at later times (Figure 17). Since, the evolution of scalar P is initially dominated 
by chemical conversion and not molecular mixing, the accuracy of the conditional dis- 
sipation model at the early stages may not be critical to the overall performance of the 
model. At the later stages, when the evolution of scalar P is mixing dominated, the 
model is adequate. For this reason, the use of equations (42 and 41) as closure model 
may yet yield reasonable results. 

2. A close comparison (of the normalized moments) of the multiscalar mapping closure 
model and the multivariate /5-pdf model (Figures 22 and 23) against inert data shows 
that the two models are quite close when the initial length scales and diffnsi vi ties 
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are similar. The mapping closure model is superior when those parameters are vastly 
different for the different scalars. 

3. The assumed /?-pdf model comes reasonably close in calculating many of the unnor- 
malized moments of the scalar joint pdf, even for the reacting case and even when 
the initial length scale and diffusivities of the participating scalars are quite disparate 
(Figures 18 - 21). It is the unnormalized moments that require closure modeling in tur- 
bulent calculations. Hence, despite the disagreement of the normalized moments under 
some circumstances, the assumed /9-pdf model appears to be adequate for engineering 
calculations. 
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Meons for Group 1 



Figure 1: Evolution of the scalar means for Group 1. pa : Cl (solid line), Cla (■ -) 
Cib (- -), Clc( ). HP : Cla (circle), Clb (square), clc(tnangle). 



Time t 


Figure 2: Evolution of the scalar means for Group 2. (j,a • C2 (solid line), C2a (- -), 
C2b (- -), C2c(— — ). up : C2a (circle), C2b (square), c2c(triangle). 
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Figure 3: Evolution. of the scalar variance for Group 1: Legend same as figure 1 


Variances for Group 2 



0.00 0.01 0.02 0.03 0.04 0.05 

Time t 

Figure 4: Evolution of the scalar variance for Group 2: Legend same as figui 






Mean scalar dissipation and chemical production rate (Group 2) 



Figure 7: Mean scalar dissipation and chemical production rate of of scalar A for Group 

2. Mean scalar dissipation : C2 (solid line), C2a (- -), C2b ( ), C2c( ). Chemical 

production : C2a (circle), C2b (square), C2c(triangle). 



Figure 8: Mean scalar dissipation and chemical production rate of scalar P for Group 1. Mean 

scalar dissipation : Cl (solid line), Cla (- -), Clb ( ), Clc( ). Chemical production : 

Cla (circle), Clb (square), Clc(triangle). 
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Figure 9: Mean scalar dissipation and chemical production rate of of scalar 1 
2. Mean scalar dissipation : C2 (solid line), C2a (- -), C2b ( ), C2c( 

production : C2a (circle), C2b (square), C2c(triangle). 


for Group 
Chemical 
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PDF of scalar A (Case 2) 
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Figure 16: DRS data vs. unmodified mapping closure model. Comparison of F(<{>a)x(^a)/^a 
for Case C2b: Data (solid line), model (circle), (a) t = 0.006, (b) t = 0.012, (c) t = 0.024 
and (d) t = 0.048. 
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Figure 17: DRS data vs. unmodified mapping closure model. Comparison of F((f)p)x{(j>p)/^p 
for Case C2b: Legend same as figure 16. 
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Comparison DRS vs. Beta— Pdf model (Case 1c) 
0.25 i 1 



Time t 

Comparison DRS vs. Beta-Pdf model (Case 1c) 



Time t 

Figure 18: DRS data vs. assumed /?-pdf model. Comparison of various unnormalized mo- 
ments for case Clc. (a) af: data (solid line), model (circle); a\ - data (- -), model (square), 
(b) (<j>i<f>- 2 ) - data (solid line), model (circle); (tffa): data (- -), mode] (square); and, (<£,$) 
- data ( ), model (triangle). 





Comparison DRS vs. Beta-Pdf model (Case 2c) 



Comparison DRS vs. Beta-Pdf model (Case 2c) 



Time t 


Figure 19: DRS data vs. assumed /i-pdf model. Comparison of various unnormalized 
ments for case C2c. Legend same as figure 18. 


mo- 
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Comparison DRS vs. Beta-Pdf model (Case 3b) 



Comparison DRS vs. Beta— Pdf model (Case 3b) 



Figure 20: DRS data vs. assumed /3-pdf model. Comparison of various unnormalized mo- 
ments for case C3b. Legend same as figure 18. 
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Comparison DRS vs. Beta-Pdf model (Case 4b) 



Comparison DRS vs. Beta— Pdf model (Case 4b) 



Figure 21: DRS data vs. assumed /9-pdf model. Comparison of various unnormalized mo- 
ments for case C4b. Legend same as figure 18. 
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Figure 22: Comparison of mapping closure model and /i-pdf model using DRS data for Case 
C3. Data (soild line), mapping-closure model ( ) and /9-pdf model (circle). 
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Figure 23: Comparison of mapping closure model and /3-pdf model using DRS data for Case 
C4. Data (soild line), mapping-closure model ( ) and /3-pdf model (circle). 
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